Vibrational spectrum of topologically disordered systems 
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The topological nature of the disorder of glasses and supercooled liquids strongly affects their high- 
frequency dynamics. In order to understand its main features, we analytically studied a simple 
topologically disordered model, where the particles oscillate around randomly distributed centers, 
interacting through a generic pair potential. We present results of a resummation of the perturbative 
expansion in the inverse particle density for the dynamic structure factor and density of states. This 
gives accurate results for the range of densities found in real systems. 
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The high frequency dynamics of glasses and super- 
cooled liquids has recently received a large amount of ex- 
perimental numerical |^-^, and theoretical 
attention. High-resolution inelastic X-ray scattering 
techniques have made accessible to experiment the region 
where the exchanged external momentum p is compara- 
ble to poj the maximum of the static structure factor. A 
number of facts have emerged from these experiments: 
a) The dynamic structure factor (DSP) S{p, uj) has a 
Brillouin-like peak for momenta up to p/po ^ 0.1 — 0.5, 
usually interpreted in terms of propagating acoustic-like 
excitations, whose velocity extrapolates to the macro- 
scopic sound velocity when p — > 0. b) The peak has 
a width P, due to disorder, which in a large variety of 
materials seems to scale as P = Ap°' , with a ^ 2 and 
A depending very slightly on the temperature, c) The 
density of states (DOS) exhibits an excess respect to the 
Debye behavior {g{ijj) <x. uP'), known as the Boson peak, 
which is most remarkable for strong glasses. 

Prom the theoretical viewpoint, the challenge is to ex- 
plain the above features of the scattering spectra. Por 
a dense system like a glass, the short time dynamics 
is naturally interpreted as a consequence of vibrations 
around a quenched disordered structure. Indeed, molec- 
ular dynamics simulations have shown that harmonic vi- 
brations are enough to describe the dynamic structure 
factor (DSP) |^ and specific heat Q|. Nevertheless, the 
analytical problem is very hard even in this first approx- 
imation, because it involves a matrix (the Hessian) with 
random elements. So the study of vibrations in glasses 
is related to the general problem of analyzing the sta- 
tistical properties of large random matrices This 
is a venerable problem, with applications that include 
the theory of nuclear spectra, conductivity in alloys, and 
many others p^ . We should distinguish between ma- 
trices obtained as a result of random perturbations of a 
reference lattice, and those for which no reference lat- 
tice can be defined |Q. Lattice-based random matrices 
arise typically in problems of the solid state (Anderson 



localization, transport in mesoscopic systems, etc.). A 
number of well-known approximations have been devel- 
oped to deal with these kind of matrices: this approxima- 
tions can be derived by considering the disordered part 
of the random matrix as a perturbation of the solvable 
crystalline problem [jl^ and partially resumming the per- 
turbative expansion. On the other hand, off-lattice ran- 
dom matrices arise in many problems, such as amorphous 
semiconductors, instantaneous normal modes in liquids, 
and the vibrational excitations of glasses. These matri- 
ces lack a natural separation between a solvable ordered 
part and a disordered perturbation. As a consequence, 
the theoretical tools for their study have started to be 
developed much more recently |15|, and are best suited 
for low densities, with emphasis on the DOS rather than 
the DSP. Purthermore, the Debye behavior of the DOS is 
not reproduced by the above methods, even in the glassy 
phase (see Cavagna et al. in ]T5[|). 

The purpose of this Letter is to address the prob- 
lem of harmonic vibrations in glasses and to find a self- 
consistent integral equation. This equation is somewhat 
analogous to the lattice-based approximations mentioned 
above, but is appropriate for the off-lattice matrices that 
arise in this case Our approach is based on an ex- 
pansion in 1/p [p is the number of particles per unit 
volume) §,§. 

Somehow complementary to the present approach is 
a recent modification of Mode Coupling Theory to 
describe excitations around a quenched structure, which 
has been applied to a hard-sphere fiuid. Also, lattice- 
based models have been proposed JlO| , p| , but it turns 
out ||l2| that they definitely miss the p^ behavior of the 
peak width. Other differences between the lattice and 
off-lattice cases regarding the DOS will be discussed. 

We consider a model of particles oscillating harmoni- 
cally around (disordered) equilibrium positions. Por sim- 
plicity we assume that they can move only along a direc- 



tion u: Xj(t) 



+ Uipj{t), J = 1,, 



. N. The vibra- 
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tional potential energy is then V{{(pi}) = \ ^ fi'^T^ ~ 
x'^'^){ipi — (fij)"^, where / is the second derivative of the 
pair potential. Of interest is the spectrum of the dynam- 
ical matrix (Hessian) 



N 
k=l 



(1) 



Neglecting quantum effects, in the harmonic approxi- 
mation the DSF can be obtained as (see e.g. |0]) 



S{p,uj) 



kBTp" 



E 



5{lo-u^), (2) 



where e„_i is the z-th component of the 71-th eigenvector 
of the Hessian (|l]) , the eigenfrequencies are the square 
root of the corresponding eigenvalues, T is the tempera- 
ture and the overline stands for the average over {x°'^} 
(hereafter we set the mass m = 1). 

The model is defined by the probability distribution 
of the {iCj'^} and by the function /(r). Obviously, in a 
real system the are highly correlated. However, as 

discussed in detail in |9|[l^, it is a good approximation 
to take {x'^^} uncorrelated while at the same time using 
an effective interaction which weights the physical / with 
the radial distribution function g{r). For a given /, the 
only free parameter is p. 

As usual, we consider the resolvent G'(p, z) = 

J2jk exp[ip • {x'j'^ - x''f^)][z - M]~^, since the DSF can be 
obtained as S{p,Lu) = -2kBTpHmG{p,uj^ + iO+)/(a'7r). 
Our aim is to compute G using the appropriate self- 
consistent equations. Although here there is no sepa- 
ration between an ordered reference state and a random 
perturbation, it has been observed 1^ that in the limit 
of infinite density the spectrum coincides with that of an 
elastic homogeneous medium, since the resolvent is sim- 
ply Gip,z) = [z- A(p)]-i (where A(p) = p[/(0) - /(p)] 
and f{p) is the Fourier transform of /(r)). A systematic 
perturbative expansion in 1/p has been constructed in 
ref. §. Technically, an expansion is found for the resol- 
vent, which is written as 



Gip,z) 



- X{p) - 5](p, z) 



(3) 



The self-energy I](p, z) — S'(p, z) + iS"(]3, z) is then ex- 
panded in powers of 1/p ||^,^ in the relevant region where 
p = 0{z). At infinite density, S = and the dispersion 
relation is (^(p) — \J A(p), which is linear at small p. The 
real self-energy E'(p, z) renormalizes the dispersion re- 
lation, while the imaginary part gives the peak width: 
r = S"(p, a;p)/aj]3, where Wp is the position of the maxi- 
mum of the DSF. 

We have reformulated the \ j p expansion [ p7[ in a way 
that reduces the number of diagrams and allows to iden- 
tify those with the simple topology of Fig. Topologi- 
cally, these diagrams are exactly those considered in the 



usual lattice CPA and in other self consistent approxi- 
mations. The sum of this infinite subset is given by the 
solution of the integral equation: 



S(p, z) 



(2^ 



P /(?)-/> -g) G{p,z), (4) 



where the resolvent is given by Eq. (|^) . The solution gives 
us the resolvent, and hence the DSF and DOS (Eq. ^ 
below) . 
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FIG. 1. The diagrams of the 1/p expansion which are taken 
into account in our approach. The numbers correspond to the 
particle-label repetitions 

We are interested in the solution of Eq. (^) for different 
values of z and p. To be definite, we consider an explicit 
case where the function /(r) has a simple form, namely 
/(r) = exp[— r^/ (2(t^)]. This is a reasonable first approx- 
imation for the effective interaction We shall take 
a as the unit of length and set po = which is a rea- 
sonable choice for po for this Gaussian /(r), as discussed 
in JlBf . In this particular case we will solve numerically 
the self-consistence equation. We will also evaluate by 
simulation (using the method of moments JlSt ) the exact 
DSF and DOS by computing the resolvent for concrete 
realizations of the dynamical matrix, considering a suf- 
ficiently high number of particles so that finite volume 
effects can be neglected. These numerical results will be 
supplemented by analytic results, that are /-independent 
and can be obtained in the limits p 00 and p — *■ 0. 

The infinite momentum limit is particularly interesting 
because of the remarkable result § that the DOS g(uj) 
can be written as 



lim 

p — ^00 



(5) 



We easily find that in this limit Eq. (||) can be written 
as {g{z) = limp^oo Gip, z) and A = (27r)-3 fp[q) d^g); 



-±- = '--f{o)-Ag{z)- 
py{z) p 



27r) 



rf\q)G{q,z) (6) 



A simple approximation consists in neglecting the last 
term in the r.h.s. of (|^), which is reasonable at large z. 
This approximation implies a DOS which is semicircular 
as a function of cj^, with width proportional to ^fp and 
centered at = pf{0)- Translational invariance also 
requires low-frequency modes. These are given by the 
neglected term, and in fact it is easy to show that at 
high density it produces a Debye spectrum which extends 
between zero frequency and the semicircular part. 

In the limit p — s- 0, the leading contribution to E" 
comes from q 3> p in Eq. (|), where G{q,z) w g{z), so 
we can write for the peak width T{p) « Tq{p), where 



2 



r,(rt.,,^/^[/,„-/(p-,f. (7) 

The integral is of order , so if the spectrum is Debye- 
Uke for small frequencies, we get T(p) ~ . 

These considerations are verified by the numerical so- 
lution of the Gaussian case, which are shown in Fig. |^ 
for p = 1.0 fj""^ together with the results for the simula- 
tions Note the good agreement, to be expected for 
high-densities, and how, for large p, S^^''{p,uj) (Fig. ||, 
top) tends to the DOS. The DOS from the self-consistent 
equation (Fig. H,bottom) also agrees very well with the 
results from simulations, and is a big improvement over 
the first term of the expansion in powers of p~^. The two 
contributions (Debye and semicircle) mentioned above 
can be clearly identified. Our approximation fails in 
reproducing the exponential decay of the DOS at high 
frequencies, which is non perturbative in 1/p [p^ and 
corresponds to localized states. 
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FIG. 2. Top: dynamic structure factor as obtained from 
Eq.(^ (full line) and from simulations |^ (dashes). Bottom: 
DOS divided by tj^ (Debye behavior) as obtained from Eq. (^) 
(full line), simulations (dashes), and first order in the 1/p 
expansion (dots). 

Next in Fig. ^ we plot the linewidth as a function of p as 
obtained from Eq. (||) . Notice that we recover the behav- 
ior predicted from the first two non-trivial terms in the 
expansion m powers of [§: the linewidth IS propor- 
tional to at small p (also predicted by the argument 
above), then there is a faster growth and finally it ap- 
proaches to a constant as S^^^ {q, uj) starts to collapse onto 
the DOS. The inset shows that the contribution Eq. (Q) 
is indeed dominant at small p. However, rigorous p'^ scal- 
ing is found only for very small momenta {p/po < 0.1), 
while experiments are done at 0.1 < p/po < 1. In this 
crossover region, our model predicts deviations from p^ , 
which are probably hard to measure experimentally. In 
any case, the effective exponent is certainly less than 4, 



in contrast with lattice models and consistent with ex- 
perimental findings. Similar conclusions can be drawn 
from mode couphng theory (Fig. 8 of 0). 
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FIG. 3. Peak width vs. p, for p — 0.6a~^. The inset shows 
that ro(p) is the dominant contribution at small p. 



Finally, let us discuss the low density case. As men- 
tioned, Eq. (H) tells us that the position of the center of 
the semicircular peak in the DOS goes to as ~ p, 
but with a width that only decreases as y/p. So at low 
enough densities the DOS develops a tail that extends 
into negative values of w^. In our Gaussian case, this 
starts happening at pc ~ O.Slcr^''. Since our original 
Hessian is positive definite, one could regard such an in- 
stability as just a failure of the 1/p expansion. However, 
the structure that Eqs. (Q) and (||) (wrongly) predict from 
the hybridization between the Debye and semicircle spec- 
tra for p ^ Pc is extremely suggestive of a Boson-Peak. 
In fact, it can be shown |l^ that for p Pc, and very 
small values of w a cross-over occurs from the Debye 
behaviour at low w to a lu^^^. This structure (Fig. |[ bot- 
tom) can be compared with the Boson peak predicted by 
the lattice integral-equations |l^,|l^ , close to their insta- 
bilities (see Fig. ||, top). The Boson peak predicted by 
the off-lattice equations is a subtle feature of the spec- 
trum that can only be seen when dividing the DOS by 
the Debye behavior. It is also striking that this off-lattice 
Boson peak shifts towards uj — when the instability is 
approached, as the experimental peak does upon heating. 
On the other hand, it is fair to say that the Boson peak of 
the lattice integral equations is just the first maximum of 
the DOS (which experimentally is not very temperature- 
dependent), shifted to low ui by the division by w^. 

In conclusion, we have presented a random-matrix ap- 
proach to study an off-lattice model for the vibrational 
dynamics of glasses, in which particles oscillate harmon- 
ically along a given direction around quenched equilib- 
rium positions. At variance with disordered systems de- 
fined as a perturbation of a reference crystal, our sys- 
tem lacks a natural separation between a solvable and 
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a random-perturbation term. Nevertheless, an expansion 
in 1/p can be obtained ||^. Here we have reformulated 
this expansion in a way that allows us to resum an infi- 
nite subset of terms, obtaining an integral self-consistent 
equation for the dynamic structure factor and the 
DOS. This equation predicts a.p^ scaling of the linewidth 
for a generic potential, as well as the appearance of a 
Boson-peak-like structure at low densities, linked to a 
mechanical instability. Since the dynamic structure fac- 
tor tends to the DOS at large momentum, we argue that 
Boson-peak-like features found experimentally in the for- 
mer at large p |^,^ are just another manifestation of the 
Boson peak. We compared the results of eq. with nu- 
merical calculations for a gaussian potential and found 
satisfactory agreement at densities comparable to those 
of real glasses. At low densities, the resummation finds 
a spectrum at negative squared frequencies, which is un- 
physical for this potential. However, for p ^ Pc & struc- 
ture appears which is definitely different from those found 
in lattice models and extremely remindful of the 

experimental Boson peak. 




FIG. 4. Top: of a;) (full line) and g(u))/uj'^ (dashes) of the 
model of Ref. |_0|, as predicted by the single-link CPA, at 
A = 1. Bottom: g{ui)/9 (full line) and g{uj)/uj^ (dashes) 
for our model at p = 0.5a~^. Note that for the values of 
the parameters shown here, the lattice model is nearer to its 
critical point ((Ac - A)/Ac ~ 0.22 while {p- pc)/pc ~ 0.61). 

These results can be extended in several ways (poten- 
tials with instabilities at low densities, non-coUinear vi- 
brations, detailed study of the correlated case); work in 
these directions is in progress. 
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